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Abstract 

Bosonic atoms trapped in an optical lattice at very low temperatures, can be modeled by the 
Bose-Hubbard model. In this paper, we propose a slave-boson approach for dealing with the Bose- 
Hubbard model, which enables us to analytically describe the physics of this model at nonzero 
temperatures. With our approach the phase diagram for this model at nonzero temperatures can 
be quantified. 



PACS numbers: 03.75.-b,67.40.-w,39.25.+k 



I. INTRODUCTION 



The physics of the Bose-Hubbard model was the subject of intensive study for some years 
after the seminal paper by Fisher et ai, which focused on the behavior of bosons in a disor- 
dered environment More recently it has been realized that the Bose-Hubbard model can 
also be applied to bosons trapped in so-called optical lattices 0, and mean- field theories 
[3, 0, and exact diagonalization jof have been succesfully applied to these systems in one, 
two and three dimensional systems. The experiments performed by Greiner et al. have 
confirmed the theoretically predicted quantum phase transition, i.e., a phase transition in- 
duced by quantum fiuctuations, between a superfluid and a Mott-insulating phase. A review 
of the work carried out in this field has been given by Zwerger Strictly speaking the 
above mentioned quantum phase transition occurs only at zero temperature j9|. At nonzero 
temperatures there is a 'classical' phase transition, i.e., a phase transition induced by thermal 
fiuctuations, between a superfiuid phase and a normal phase and there is only a crossover 
between the normal phase and a Mott insulator. It is important to mention here that a 
Mott insulator is by definition incompressible. In principle there exists, therefore, no Mott 
insulator for any nonzero temperature where we always have a nonvanishing compressibil- 
lity. Nevertheless, there is a region in the phase diagram where the compressibillity is very 
close to zero and it is therefore justified to call this region for all practical purposes a Mott 
insulator Qualitatively this phase diagram is sketched in Fig. ^for a fixed density. This 
figure shows how at a sufficiently small but nonzero temperature we start with a superfiuid 
for small positive on-site interaction U , we encounter a phase transition to a normal phase 
as the interaction strength increases, and ultimately crossover to a Mott insulator for even 
higher values of the interaction strength. We can also incorporate this nonzero temperature 
behaviour into the phase diagram in Fig. |21 This figure shows how at zero temperature we 
only have a superfiuid and a Mott insulator phase, but as the temperature is increased a 
normal phase appears in between these two phases. 

The aim of this paper is to extend the mean-field approach for the Bose-Hubbard model 
to include nonzero temperature effects and make the qualitative phase diagrams in Figs. ^ 
and 121 more quantitative. To do that we make use of auxiliary particles that are known as 
slave bosons. The idea behind this is that if we consider a single lattice site, the occupation 
number on that site can be any integer. With each different occupation number we identify 
a new particle. Although this means that we introduce a lot of different new particles, the 
advantage of this procedure is that it allows us to transform the on-site repulsion into an 
energy contribution that is quadratic in terms of the new particles. Because we want to be 
able to uniquely label each different state of the system, the new particles cannot indepen- 
dently be present at each lattice site. That is why we have to introduce a constraint. Using 
this we derive within a functional-integral formalism an effective action for the superfluid 
order parameter which depends on the temperature. The equivalence with previous work at 
zero temperature is demonstrated. 

The outline of the paper is as follows. In Sec. Owe introduce the slave-boson formalism 
and derive an effective action for the superfluid order parameter. In Sec. Illll we present the 
zero and nonzero temperature mean-fleld results. The remainder of the paper is devoted to 
the effect that the creation of quasiparticle-quasihole pairs have on the system. 
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II. SLAVE-BOSON THEORY FOR THE BOSE-HUBBARD MODEL 



In this section we formalize the above introduced idea of the slave bosons. We rewrite 
the Bose-Hubbard model in terms of these slave bosons within a path-integral formulation 
and derive an effective action for the superfluid order parameter, which then describes all 
the physics of our Bose gas in the optical lattice. 

The slave-boson technique was introduced by Kotliar and Ruckenstein [lO], who used it 
to deal with the fermionic Hubbard model. A functional integral approach to the problem 
of hard-core bosons hopping on a lattice has been previously put forward by Ziegler [11] and 
Fresard Let us first shed some light on this slave-boson formalism. We consider a single 
site of our lattice. If the creation and anihilation operators for the bosons are denoted by a| 
and a- respectively, we can form the number operator iVj = ajcij, which counts the number 
of bosons at the site i. In the slave-boson formalism, for any occupation number a pair 
of bosonic creation and annihilation operators is introduced that create and annihilate the 
state with precisely that given integer number of particles. The original occupation number 
states I Hi) are now decomposed as | n^^n], . . .), where nf is the eigenvalue of the number 
operator nf = (a")"'"a" formed by the pair of creation {a'^)^ and annihilation af operators 
that create and annihilate bosons of type a at the site i. As it stands, this decomposition is 
certainly not unique. For example, the original state 1 2) could be written as | 0, 0, 1,0,.. .) or 
as 1 0, 2, 0, . . .). Our Hilbert space thus greatly increases. To make sure that every occupation 
occurs only once we have to introduce an additional constraint, namely 



fl° = 1 (1) 

for every site j. This constraint thus makes sure that there is always just one slave boson 
per site. Because in the positive U Bose-Hubbard model bosons on the same site repel each 
other, high on-site occupation numbers are disfavored. It is therefore conceivable that a 
good approximation of the physics of the Bose-Hubbard model is obtained by allowing a 
relatively small maximum number, e.g. two or three or four, of bosons per site. 
As is well known, the Hamiltonian of the Bose-Hubbard model reads, 

H = -Y^ a^Afij «i • (2) 

Here denotes the sum over nearest neighbours, are the hopping parameters, and n 
is the chemical potential. Using our slave-boson operators we now rewrite Eq. Q into the 
form 

(iJ) a, 13 i a 

i a 

with the additional constraint given in Eq. (^J. We see that the quartic term in the 
original Bose-Hubbard Hamiltonian has been replaced by one that is quadratic in the slave- 
boson creation and annihilation operators, which is the most important motivation for the 
introduction of slave bosons. 
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Now that we have introduced the slave-boson method and derived its representation of 
the Bose-Hubbard model, we want to turn the Hamiltonian into an action for the imaginary 
time evolution. Using the standard recipe 0, 3| we find 



"^0 I i a/3 i 



-EE V^^/(TTliarratU,a^^\a^r \ , (4) 

where M is a diagonal matrix that has as the a*'* diagonal entry the term hd/dr — a// + 
a{a — l)f//2, and (3 = l/ksT is the inverse thermal energy. The real valued constraint field 
A enters the action through, 

l[6{J2n? - 1) = / d[X]e^^^o'T:.^dr){T:.nT-iyr_ 

i a 

Although we have simplified the interaction term, the hopping term has become more 
complicated. By performing a Hubbard-Stratonovich transformation on the above action 
we can, however, decouple the hopping term in a similar manner as in Ref. 0. This 
introduces a field $ into the action which, as we will see, may be identified with the superfluid 
order parameter. The Hubbard-Stratonovich transformation basically consists of adding a 
complete square to the action, i.e., adding 

/''^^E (^"^i - E v^^(«r')x j (^<^'.-E^^<^'«)*) • 

Since a complete square can be added to the action without changing the physics we see that 
this procedure allows us to decouple the hopping term. We also perform a Fourier transform 
on all fields by means of a°'(r) = (l/i/iVT^) J2kn'^i<.n^^^^^'~'^"^^ ■ If "^^ also carry out the 
remaining integrals and sums we find 



k,n ^ k,q n,n' 

+EKJ*^"'(-o<n- E tIt^ I f E v^k:.^' 

k,n k,k' ,n,n' V^V.A/^ \^ ^ 

V a 



n+n')*«k,n "^k'^ 



n+n' 

(6) 

where the matrix M{iujn) is related to the matrix M in Eq. (0} through a Fourier transform. 
Furthermore, A = (Ao,o/V^s^/3), Ck = '^t ■ cos [kj a), where a is the lattice constant of 
the square lattice with Ng lattice sites. For completeness we point out that the integration 
measure has become 

I d[{ay]d[a-] = I n^[«n)1rfKn];i- (7) 

k,n ^ 
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In principle Eq. ® is still an exact rewriting of the Bose-Hubbard model. As a first 
approximation we soften the constraint by replacing the general constraint field Ai(r) with 
a time and position independent field A. By neglecting the position dependence we enforce 
the constraint only on the sum of all lattice sites. Doing this we are only left with the Ao,o 
contribution in Eq. (jH)), which can then be added to the matrix M. The path-integral over 
the constraint field reduces to an ordinary integral. So we have, 

5[$*,$,(a°)*,a",A] = ^0 + 5/ (8a) 
where, 

^0 = ^iV,n/?A + 5^5^|ek|$k,nr + (a^,J* M^^^^iUn) <„} = So^"" + ek|$k,n|', (8b) 

a,f3 k,n k,n 

The matrix M'^^iiun) = dap^—ihjJn — iX — a/i + a(a — l)f//2), and 

(8c) 

The crucial idea of Landau theory is that near a critical point the quantity of most 
interest is the order parameter. In our theory the superfluid field $ plays the role of the 
order parameter. Only $0,0 can have a nonvanishing expectation value in our case and, 
therefore, we can write the ground-state energy as an expansion in powers of $0,0 > 

Eg{<l>o^o) = ao{a, U, fi) + 02(0, U, fi)\<l>o,o\^ + O(|$o,o|'), (9) 

and minimize it as a function of the superfluid order parameter $0,0 • We thus find that 
($0,0) = when 02(0;, U, ^) > and that ($0,0) 7^ when 02(0, U, fi) < 0. This means that 
02(0;, f/, /i) = signals the boundary between the superfluid and the insulator phases at zero 
temperature and the boundary between the superfluid and the normal phases at nonzero 
temperature. Therefore we are going to calculate the effective action of our theory up to 
second order in The zeroth-order term in the expansion of the action in powers of the 
order parameter gives us the zeroth-order contribution Qq to the thermodynamic potential 
Q. We have, 



g-/3f2o 



/ n (n^[«n)1^Kn]^) e-^o-M (10) 



a \k,ri 

From this it follows that 



-pno = -tN,(3\ + iV. ^ log (1 - e-^^^""(°)) , (11) 

and M"°(0) = (— iA — a/i + a{a — 1)U /2). Next we must calculate {S'j) where (■ ■ ■ ) denotes 
averaging with respect to 5*0, i.e., 

(^) = (n^[Kn)1rfKn]^) A[(a")*,a1e-^o-M (12) 
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Once we have this contribution, we automatically also find the dispersion relations for the 
quasiparticles in our system as we will see shortly. For small $ we are allowed to expand 
the exponent in the integrand of the functional integral for the partition function as 

^-sm ^ ^-iSo+s,)/n ^ g-5o/a^i _ Sj/h+^iSi/hf). (13) 

It can be shown that the expectation value of Sj vanishes. The second order contribution 
is found to be, 

<^') = 2 E 4]|^E(" + l)<K:k',n+n')*<:k',n+n')(Kn)*<n)- (14) 
k,k', n,n' a 

One of the sums over the Matsubara frequencies a;„ can be performed and the sum over 
k' produces an overal factor A^^. We thus find 

4^E(" + i)-^^ — - — ^' (15) 

" hd -thcon - H + all ^ ' 

k,n a 

where we have defined the occupation numbers = {{af)*af) that equal 

= -^^ . (16) 

exp |/3 — an + ^a{a — l)U) } — 1 

Having performed the integrals over the slave-boson fields to second order, we can exponen- 
tiate the result to obtain the effective action for the order parameter 

= ^W-^E$k,nG'-'(k,^a;„)$k,nj , (17) 
where we have defined the Green's function 

-hG-\k,^u;:,)= (ei.-e^^(a + l) J^" " |. (18) 

\ —inuJn — Ai + OiU J 

This result is one of the key results of this paper, which is correct in the limit of small 
$k,n- If we want to make the connection with the Landau theory again, we can identify the 
a2{(y,U,n) in Eq. with G~^(0,0)//?. In Sec. Illll we analyse this further. 



A. Mott insulator 

In the Mott insulator where hq = |($o.o)P = 0, the thermodynamic potential is now 
easily calculated by integrating out the superfiuid field. In detail 



= J d\ exp < —P^Iq — lo; 



-ihWn — + C(U 



..(19) 
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At this point we perform a saddle point approximation for the constraint field A. This 
implies that we only take into account that value of A that maximizes the canonical partition 
function. If we now thus minimize the free energy with respect to the chemical potential 
and the constraint field, we get two equations that need to be solved. The first is dfl/dX = 
and reads, 

iV. f 1 - E - ^ E ^(k^ ^-n) = 0. (20a) 

\ a / ' k,n 

In a mean-field approximation the last term is neglected, and this equation tells us that the 
sum of the average slave-boson occupation numbers must be equal to one. This reflects the 
constraint of one slave boson per site. The second equation follows from —dQ/dfi = N and 
gives 

iV. E + ^ E ^(k, ^u:n)^^^^p^ = N. (20b) 

a ^ k,n ^ 

This equation shows how the particle density can be seen as the sum of terms an" and a 
correction coming from the propagator of the superfluid order parameter. The latter is again 
neglected in the mean-field approximation. 



B. Superfluid phase 

In the superfiuid phase the order parameter |$o,oP has a nonzero expectation value. We 
find this expectation value by calculating the minimum of the classical part of the action, i.e., 
— ^G~^(0, 0)|$o,oP + '24|^o,o|^- This minimum becomes nonzero when —hG^^{0,0) becomes 
negative, and is then equal to 

I(*o,o)r = ^^^no (21) 

In appendix 1X1 we calculate the coefficient of the fourth order term |$o,o|^- We ap- 
proximate the prefactor to the fourth order term, which in general depends on momenta and 
Matsubara frequencies, with the zero-momentum and zero-frequency value of so that the 
approximate action to fourth order becomes, 

S = hpno-h^<^l^^fi-\k,iuJn)^i,^n + ai E E '^k,n'^k',n''^k",n"$k+k'-k",n+n'-n" (22) 
k,n k,k',k" ■n,n',n" 

We now write the order parameter as the sum of its expectation value plus fiuctuations, i.e., 
"^0,0 ^ y/^o ■ y/Nsh(3 + $0,0 and a similar expression for $o o- P^^ ^^is into the action 

and only keep the terms up to second order, the contribution of the fourth-order term is 
given by 

a^no E ('^k.n'f -k,-n + ^K,n^Kn + '^k,n'^'-k,-n) • 
k,n 
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There is also a contribution —hG ^(0,0) no from the second-order term. To summarize, in 
the superfluid phase we can write the action Eq. to second order as 



k,n ' ' 

^ 2^a4^o -G-\-k,-iuj,,)+4ha4noJ ■ ^^"^^ 

Integrating out the field $k,n we find the Bogoliubov expression for the thermodynamic 
potential in the superfiuid phase, 

= J dXexp{-pno + noG-\O,0) -TT[\og{-hpG~^)]} (24) 



III. MEAN-FIELD THEORY 



In this section, we apply the theory we have developed in the previous section. First, 
using the Landau procedure, we reproduce the mean-field zero-temperature phase diagram. 
We then study the phase diagram at nonzero temperatures. To do so we calculate the 
compressibillity of our system as a function of temperature, showing how for fixed on-site 
repulsion U the Mott insulating region gets smaller. By also looking at the condensate 
density as a function of temperature, we get a quantitative picture of what happens at fixed 
on-site repulsion U. The nice feature is that all our expressions are analytic. Next, we 
consider our system at zero temperature again and we study at the mean-field level the 
behaviour of the compressibillity as we go from the superfiuid phase to the Mott insulating 
phase. What we find is consistent with the general idea that the quantum phase transition 
between the Mott insulator and the superfiuid phases belongs to different universality classes 
depending on how you walk through the phase diagram (cf. Ref. [^). We then obtain an 
analytic expression for the critical temperature of the superfiuid-normal phase transition in 
the approximation of three slave bosons, i.e., up to doubly-occupied sites. Numerically we 
extend this study to include a fourth slave boson and find only slight changes to Tc. From 
the propagator of the superfiuid field we extract the dispersion relations of the quasiparticle- 
quasihole pairs and their temperature dependence. 



A. Zero-temperature phase-diagram 

From the zeros of 6*^^(0, 0) in Eq. (fTHj) . we obtain the mean-field phase diagram in the 
(yU, f/) plane. For a Mott insulating state with integer filling factor a' we have = 6a,a'- 
When this is substituted into the equation G~^(0,0) = we can find the f/(/i) curve that 
solves that equation and thus determines the size of this Mott insulating state. For given 
filling factor a' we also define Uc as the minimal U that solves the equation. Within the 
Mott insulating phase we have a zero compressibility k = dn/dfi, where n = n{fi, U) is the 
total density as determined from the thermodynamic potential. Straightforward calculation 
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gives that we are in a Mott insulating phase whenever /i lies between /i^' and /i^' where, 



fii = \ {U{2a' - 1) - 1) ± \^m-2U{2a' + l) + l. (25) 

Here we have introduced the dimensionless chemical potential p, ^ fi/zt and on-site repulsion 
strength tJ = U/zt. When Ji does not lie between any /i" and /i" the 'superfluid' density 
|('^'o,o)P will no longer be zero and the Mott insulating phase has disappeared. We have 




Oiir slave-boson approach reproduces 
For nonzero temperatures the 



drawn the zero temperature phase diagram in Fi 
here the results of previous mean-field studies [1 
equation G~^{Q, 0) = no longer describes a quantum phase transition between a superfluid 
and a Mott insulator but it describes a thermal phase transition between a superfluid and 
a normal phase. We will look into this in more detail in Sec. IIIIFI 



B. Compressibillity 

To see what happens to the Mott insulator as we move away from zero temperature we 
must look at the compressibillity as a function of temperature. Numerically we have solved 
Eq. (plj) . which gives us the occupation numbers of the slave bosons as depicted in Fig. HI 
With that we can determine the total density in the phase where the order parameter is zero. 
It is clear that within a mean-field approximation the compressibillity at zero temperature 
is exactly zero. In Fig. HI we have plotted the total density as a function of temperature. As 
the temperature is raised we see that the compressibillity, which is the slope of the curve, 
becomes nonzero. This shows that there is no longer a Mott insulator present. We also 
see that even though the slope is no longer zero it is very small indeed for low enough 
temperatures. We determine the crossover line by requiring that A{T)/kBT is of order one, 
where A(T) is defined as the difference of the quasiparticle and quasihole dispersions at 
k = 0. 



C. Superfluid density 

In a mean-field approximation the superfluid density is extracted from the action by 
finding the |($o,o)P that minimizes the fourth-order action in Eq. ()22j) . 

whenever /x is not between /i"' and and zero otherwise. We have plotted this expectation 
value in Fig. Elfor a' = 1. In this figure we see how the superfluid density grows as a function 
of n moving away from the Mott insulator phase. Our expansion of the Landau free energy 
is only valid around the edge of the Mott lobes and therefore breaks down when we go too 
far away from the Mott insulator. This can be seen in the figure as the decrease of the 
superfluid density when /i approaches and/or U. It can also be seen from the propagator 
of the superfluid field, which has poles when fi = all. For U not too far away from the 
insulating phase the figure quantitatively agrees with the ones calculated by other authors 
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D. Bogoliubov dispersion relation 



We now demonstrate that the dispersion hwk is hnear in k in the superfiuid phase and 
that the spectrum is gapless. In the superfiuid phase we can expand around the expectation 
value no = hG~^{0, 0)/2a4 of the order parameter. Up to quadratic-order this gives, 

S = hi3no - hJ2K,nG~\K ^^n)*k,n + ^4^0 ($k,n*-k,-n + ^K,nK,n + *k,n*-k,-n) ■ 
k,n k,n 

(27) 

From this we find the dispersion-relation Jtw^ in the superfiuid in the usual way. We perform 
an analytic continuation G~^(k, ia;„) — > G~^(k, 0;^) and find 

^k = h^{G-^{k,u;y,)/2 - G-^{O,0)f - {G-'{0,0)/2f. (28) 
Note that (k, a'k) = (0, 0) is a solution. Expanding around this solution in k now gives, 

fiWk^a ^^|k|, (29) 

v2 

where a is again the lattice constant. 



E. Near the edges of the Mott lobe 

If we substitute the vacuum expectation value of the order parameter back into our 
effective action, we see that the zeroth-order contribution to the thermodynamic potential 
in the superfiuid phase in mean-field approximation is given by, 

fi/?n=fi^n„-'«g-:'°-°»'. (30) 

2 04 

Prom this the particle density can be obtained by making use of the thermodynamic identity 
N — —dVt/dii. We can calculate this at T = and take the hmit // — > to show that 

the derivative of the density with respect to /x, i.c, dn/dji shows a kink for all U 7^ Uc- 
This means that only if we walk through the tip of the Mott lobes there is not a kink in 
the compressibility. In fact it's not hard to see why this is true. At zero temperature the 
roots of —hG~^{0,Qi) are by definition This means that we can write — ^^^"^(O, 0) = 
C{fi — fi'^){fi — /i"). The proportionality constant can be shown to be equal to C = 
60/ {{a'U — ii){{a' — 1)U — /i)). This then shows that the thermodynamic potential is, 

4 04 

Remembering that the density is the derivative of the thermodynamic potential we see that 
the second derivative of the thermodynamic potential with respect to ^ can show a nonzero 
value upon approaching the Mott lobe. Since in the Mott isolator the density is constant 
and equal to a' we have shown the existence of a kink in the slope of the density for all paths 
not going through the tip of the Mott lobe. This causes the difference in the universality 
class of the quantum phase transition. 
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F. The superfluid-normal phase transition 



In this subsection, we show that it is possible to obtain an analytical expression for the 
critical temperature of the transition between superfluid and normal phases as a function 
of U, for values of U below the critical U of the zero-temperature superfluid-Mott insulator 
transition. The analytical result is obtained if we include occupations up to two per site, 
i.e., three slave bosons or occupation numbers n^, n^, n^. Along similar lines can be found 
numerically if more slave bosons are included. We have carried out this procedure for the 
case of adding a fourth boson (triple occupancy) and find only modest quantitative changes. 

If we restrict the system to occupancies 0, 1 and 2, and fix the total density n = N/Ns 
at 1, the occupation numbers n^, and should obey the following relations if we neglect 
fluctuation corrections (cf. Eq. ^B^ ): 

n° + + = 1, (32) 

and 

+ 2n^ = 1. (33) 

The n" are furthermore given by Eq. enabling us to eliminate A and express and 

in terms of n^. We obtain 

(ni + l)exp(/?/i) -ni' ^ ^ 

and 

(ni + l)exp(/?(f/-/i))-ni' ^ ^ 

The constraints in Eqs. (jH^ and (jHHI) immediately lead to n'^ = n^, so that, according to 
Eqs. flMj) and we must have fi = U/2. We notice that at this level of approximation, 
we obtain a slight discrepancy with the result from Sec. IIII Al that at zero temperature the 
critical value of U of the superfluid-Mott insulator transition, which is the limiting U for 
the superfluid-normal transition that is addressed here, is according to Eq. (j25j) with a' = 1 
determined hj p, = {U — l)/2 (15>] . 

As argued above the criticality condition for the superfluid-normal transition is obtained 
by putting 6*^^(0,0) = 0. Restricting the sum in the right-hand side of Eq. (frR|) to a = 
and a = 1, we obtain llf 



1 = {n' -n')+- {n' - n°) . (36) 

Since the relation between /i and U is fixed by Eqs. (j32|) and (j33p . and nP and can be 
expressed in as n° = = (1 — n^)/2, the criticality condition Eq. (|36p results in a 
remarkably simple relation between and U at T^, namely = {(J + 3)/9. Using this in 
Eq. ()34j) leads to the following analytic formula for Tc = Tc/zt for the superfluid-normal 
transition: 



(t/-24)([/ + 3) 
(f/-6)(f/ + 12) 



(37) 
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It is straightforward to generalize this procedure to arbitrary integer density ol while allowing 
occupation numbers n" ~^ ^rf , only. The result is 



k^f:' = I loj 



{U -8{2a' + 1)){U + {2a' + 1)) 
_{U - 2{2a' + l))(f7 + 4(2a' + 1)) 



(38) 



The critical temperature T^. for integer filling factor n = N/Ng = 1, i.e., Eq. (jHTj) . is 
plotted in Fig. IHl The overall qualitative behavior is as one would expect (cf. Fig. [l]). A 
few finer details appear to be less satisfactory. For instance, vanishes for If = Q, whereas 
we would expect this to coincide with the mean-field result for Uc for the superfluid-Mott 
insulator transition for the first Mott lobe, i.e., Uc = 5.83 obtained from Eq. ()25|1 with 
a' = 1. We note that the discrepancy is not large and is even smaller for the higher Mott 
lobes. Indeed U{Tc ^ 0) = 2(2a' + 1) versus Uc = {2a' + 1) + ^/i2a^TW^. Another 
feature is the maximum in the Tc{U) curve (cf. Fig. ^and 3]). Both features mentioned are 
caused by the fact that the two conditions Eqs. and (jHHj) are strictly enforced, whereas 
they become less appropriate for small U. The exact solution for four slave bosons on 
a four site lattice for small U shows that a better result may be obtained if a fourth boson 
occupation number is included in our approach. The set of equations to be solved then 
becomes, again for n = 1, 

n° + + + = 1 (39) 
+ 2n^ + 3r2^ = 1 (40) 
3 (^3 _ ^2^^ ^ ^^2 _ ^1) + 1 (^1 _ ^ ^ _ (41) 



fi — 2U ^ ' /i — f/^ ' n 

Again n'',n^, and can easily be expressed in terms of n^, but no exact solution appears 
to be possible in this case. However, we have managed to find solutions numerically. The 
results for are depicted in Fig. IHl and show fairly little quantitative change compared to 
the analytical result Eq. ()37p. In particular, Tc still vanishes for U ^ 6, and the maximum 
is still there, although shifted to a lower U ^ 1.8 compared to U = 2.15 for Eq. (|H7|). It is 
satisfactory to find that for the higher values of U, starts to increase rapidly towards 1, 
signalling the approach of the Mott-insulator phase, whereas is almost negligible (< 1%) 
already for U ^ 3, supporting a description in terms of 3 slave bosons only ^18i] . 



G. Quasiparticle-quasihole dispersion relations 

Consider now the propagator G^^(k, cu), given by 



At zero temperature and for a given integer filling factor a', we have in a mean-field approx- 
imation that n° = 6o,n' and we retrieve the previously found result for the quasiparticle- 
quasihole dispersions In this case the real solutions of hw follow from a quadratic equation 
G^^(k,ujn) = 0. At nonzero temperature the occupation numbers in general are all nonzero 
and there will be more than just two solutions for hw. In the set of solutions there are still 
two solutions that correspond to the original single quasiparticle and quasihole dispersions. 
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The physical interpretation of the other solutions is that they correspond to the excitation 
of a higher number of quasiparticles and quasiholes. In Fig. [7| we show the three low lying 
excitation energies for k = at a temperature of zt(3 = 10. To obtain an analytic expression 
for the single quasiparticle-quasihole dispersion we only take into account the two terms in 
the sum in Eq.(|18|) which have numerators ~^ — and n° — n° These correspond 
to processes where the occupation of a site changes between a' — l,a' and a' + 1. We find 

''''^ = -^ + ^ + ^e^{a'n^'"' - n"' + {a' + l)n"'+') 

±-Ju^ + 2(a'n"'-i - (1 + 2a')n^' + (1 + a')n'^'+^)Ue]^ + + - (1 + a')n^'+^Ytl. 

(43) 

In Fig. ISlwe have plotted these dispersions at k = as a function of U. Comparison with 
Fig. m shows that Eq. (jl3|l gives an appropriate description of the single quasiparticle- 
quasihole dispersions. As can be seen from Fig. |Slthe tip of the lobe moves to smaller U 
as a function of increasing temperature. This can be understood because that point now 
describes the superfiuid-normal phase transition (cf. Figs. [H IHI). In Fig. El we show how 
the superfiuid-normal boundary in the Ji — U plane evolves for nonzero temperatures. If we 
define the gap as the difference between the two solutions at k = 0, we find that the gap 
grows bigger as the temperature increases. As we have seen in Sec. IIIIBI it is incorrect, 
however, to conclude from this that the region of the Mott insulating phase in the /i-f/ phase 
diagram grows as temperature increases. As mentioned previously, strictly speaking there 
is no Mott insulator away from zero temperature and at nonzero temperatures there is only 
a crossover between a phase which has a very small compressibillity and the normal phase. 



IV. FLUCTUATIONS 



In this section we make a first step towards the study of fluctuation effects and derive 
an identity between the atomic Green's function and the superfiuid Green's function in Eq. 
p8|) . This we then use to calculate the atomic particle density. In appendix ^ we show 
that the easiest way to calculate the density is by making use of currents that couple to the 
atomic fields. We start with the action of the Bose-Hubbard model 



S[a* , o] = f dr a* | /i— /i j — tija*aj H a*c 



(44) 



We are interested in calculating the {a*ai) correlation function. Therefore we add currents 
J* , J that couple to the a* and a fields as 

Z[J*,J] = / d[a*]d[a]exp < —So/h+ — / dT\^a*tijaj + / (ir [J*aj + a* Jj] >. 

(45) 

Here 5*0 = 5*0 [a*, a] denotes the action for tij = 0. The most important step in the remainder 
of the calculation is to a perform again a Hubbard-Stratonovich transformation by adding 
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a complete square to the action. The latter can be written as, 



Straightforward algebra yields 



L k,n ^ 

(47) 

Differentiating twice with respect to the currents gives then the relation 



1 52 

-z[r,J] 



Z[0,0]SJlJJ^, 



= «,n«k,n) = ($k,n*k,n) " ^- (48) 

J*,J=0 



This is very useful indeed since the correlator ($^^$k_n,) = —G{k, iun)- At zero temperature 
the retarded Green's function can be written as 

where the wavefunction renormalization factor is 

^ _ ^(1 + 2a') -e^+ ^U^ - 2Ue^il + 2a') + ej 
" 2^U^-2Ue^{l + 2a') + ei 

and 



e 



qp,qh 
k 



^ + |(2«' - 1) - y ± i\/4 - (4«' + 2)f/ek + (49c) 



Note that is always positive and in the limit where ?7 ^ oo we have that (1+a'). 
The quasiparticle dispersion e^^ is always greater than or equal to zero and e'^ is always 
smaller than or equal to zero. Because of this only the quasiholes give a contribution to the 
total density at zero temperature. The density can be calculated from, 

k,n k,n ^ k 

If we expand the square- root denominator of Z for small k we see that it behaves as 1/k, 
therefore in two and three dimensions we expect the integration over k to converge. In 
Fig. Uniwe have plotted the density for a' = 1 as given by the equation above. We see that 
the density quickly converges to one, but near the tip of the Mott lobe in all dimensions 
it deviates significantly from one. This result is somewhat unexpected and may be due to 
the break-down of the gaussian approximation near the quantum phase transition. A more 
detailed study of the fluctuations is beyond the scope of the present paper and is therefore 
left to future work. 
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V. CONCLUSIONS 



In summary, we have applied the slave-boson formalism to the Bose-Hubbard model, 
which enabled us to analytically describe the physics of this model at nonzero tempera- 
tures. We have reproduced the known zero-temperature results and we have computed the 
critical temperature for the superfluid-normal phase transition. The crossover from a Mott 
insulator to a normal phase has also been quantified. We have shown how thermal fluctua- 
tions introduce additional dispersion modes associated with paired quasiparticles-quasiholes 
propagating through the system. We have also considered density fluctuations induced by 
the creation of quasiparticle-quasihole pairs. These fluctuations do not average out to zero 
in the gaussian approximation. 
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APPENDIX A: HIGHER-ORDER TERMS 

If we also want to calculate quantities like the superfiuid density, we have to calculate 
the effective action up to fourth order. One way to do this is by going to higher order in the 
interaction part. Here we follow a slightly different strategy. Because we are only interested 
in the mean-field theory, it suffices to just consider $o,o terms. The effective action for $o,o 
is found from 

/■d|*; jd|*„,„l /■ d|Kt,,M<_„] /i \ 

^=y — w — yn — — '^^[-k')' '^i' 

a,k,ra 

where from Eq. (jHl) we have 

S = iNM^ + eo|$o,o|' + E«n)*^"''<n- (A2) 

af} k,n 

Note, however, that now the matrix M is only blockdiagonal and it contains off-diagonal 
terms proportional to $o,o- When we take the determinant of that matrix, you get automat- 
ically all powers in $o,o- This can be made more explicit by looking at the block-structure 
of the matrix which is 



M 



(Bo \ 

B2 

V •••/ 



(A3a) 



where 



^"=Ut?:*. . (A3b) 

with Xa = —ihujn — — afj, + a{a — l)U/2. The slave bosons can be integrated out with 
the result 

Z = l «:^fl5^exp |-1 {^NMX + 6o|$o,or) } exp | - log [det/5M]|. (A4) 
The determinant can be calculated up to fourth order in $o,o as 



(n^ \(^^\^ \cf. |2 (« + 1) , V- 4 1^ ,4(« + l)(/3 + l) 

det/3M= ^ll(3x.j (^1 + L]^l^o,o| + L L JWJ^.l^o,^ 



\a-f3\>2 (^'^(^y ' XaXa+lXf3Xf3+l 



(A5) 
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For small $0,0 we can expand the logarithm in Eq. ()A4|) by using the Taylor expansion 



log {1 - ax^ + 7x^} = ~ax^ + l/A{-2a^ + 47)^^ + 0{x^). 



Combining the latter equation with Eq. ()A4|) . we also recover that the second-order term 
in the effective action for $0,0 is given by 



eo 

k,n a 



2 



(A6) 



We determine the effective action to fourth order in the case of the first four slave bosons. 
Using the above we can readily verify that 



^0 \2 / 3 2 1 \ 1^ |2 , / ^0 n4 3 1^ ,4 



Nhfj \X3X2 X2X1 XiXoJ NshfJ X0X1X2X3 



From this we find that 04 in the case of four slave bosons is given by 
eo V^l fJ^(a + l)Y 12 



(A7) 



4 \y/NshPj \ l^XaXa+l/ XoXlX2X3 



or explicitly, 

eo 



a4 



- ^ ^^^^ - V - ^ (^^^° - (^^W^^^ 

4 1,, u 4t/ 1 12 3 

+ ^ — n (1 — n — TTT^n — - — = — — = —n 

{U-fi)fi ^ ' {U-Jifji^ {W - 2ji){2U - Jif 

12 5, 12U , 12 1 



{2U^ -3Ufi + p,^) ' ' (2f/2 - 3f//i + /i2)2 (3f/-2/i)(f/-^)2 



(A9) 



Note that in the zero-temperature limit for the first Mott lobe, when the slave-boson occu- 
pation numbers are proportional to a Kronecker delta, this result coincides exactly with the 
one previously derived in standard perturbation theory (cf. Ref. Q ). 
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APPENDIX B: DENSITY CALCULATIONS 



In this section we demonstrate for the noninteracting case the equivalence of the calcu- 
lation of the total particle density through the thermodynamic relation = —dn/dfi and 
through the use of source currents that couple to the atomic fields. We consider a system 
of noninteracting bosons described by creation and annihilation fields a*(r) and ai{r) on a 
lattice. First we calculate the generating functional Z[J*, J] for this system, 



Z[J*, ~ J d[a*]d[a] exp < ——So[a*, a] + — J dr a*tijaj + J dr {JiCii + alJi 

\ ij i 



(Bl) 



In this equation 5*0 is the on-site action, which in frequency-momentum representation typ- 
ically looks like 

So[a*, a] = E (-«^n - ak,n- (B2) 

k,n 

The hopping term is decoupled by means of a Hubbard-Stratonovich transformation, i.e., 
we add the following complete square to the action, 

E ( < - + ^E^^'-^' ) U - + ^E^^'"^^" ) ■ 
ij \ j' J \ j" / 

The atomic fields a*, a can now be integrated out. Going through the straightforward 
algebra one arrives at the following expression for the generating functional. 



Z[r, J] = J d[^*]dm exp I K,nG-'{k, Z^„)$k,n + J;„$k,n + Jk,n$k,n ' ^ 



k,n 



(B3) 



where — ^G~^(k, ia;„) = ek — {—ihwn — A*) ^ The total density may be calculated from 
this expression by first calculating the correlator (Okn^k.n) through functional differentiation 
with respect to the source-currents J, and then to sum over all momenta and Matsubara 
frequencies. We have for the first step 



—ihuJn ~ A* ~ Ck 
J*,J=0 n P k 



' (B4) 



We sec that there is a pole here at = — Ck — jJL. The density now can be calculated from 
n = (l/NshP) X]kn('^kn'^k,n)- This is thc expected result. 

On the other hand, we can also calculate the density from the thermodynamic potential 
Q, by using the relation N — —dQ/dfj, where N is the total number of particles. Doing that 
for this case we use that 

^ = ^ E IPH^n - 1^)] + log [-h(3G-\k, tUJr^)] } (B5) 

^ k,n 
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and obtain 

1 dfl 1 ^ ( h h ek 

Ng da Nghp ^ \ —ihjJn — U —ihujn — M — Ck —ihujn — A* 

k,n 

When doing the sum over Matsubara frequencies the pole at ihiUn = — A* in the first term 
the right-hand side is canceled by the second term and only the other pole at ihuJn = — Ck — 
gives a contribution. This shows the equivalence of both methods. 



(Bi 
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FIG. 1: Qualitative phase diagram for a fixed and integer filling fraction in terms of the temperature 
T and the dimensionlcss coupling constant U / zt, with superfluid (SF), normal and Mott insulating 
phases (MI). Only at T = a true Mott insulator exists. 




FIG. 2: Qualitative phase diagram in terms of the chemical potential zt and the dimensionless 
coupling constant U/zt. For nonzero temperatures a normal phase appears. 
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FIG. 3: Phase diagram of the Bose-Hubbard Hamiltonian as obtained from the mean-field zero- 
temperature hmit in the slave-boson formalism. It shows the superfluid (SF) phase and the Mott 
insulator regions with different integer filling factors here denoted by a' . The vertical axis shows 
the dimensionless chemical potential /2 = ^/zt and the horizontal axis shows the dimensionless 
interaction strength U = U jzt. 
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FIG. 5: Supcrfluid density |<&o,oP a-s a function of for various temperatures and for U/zt = 10. 
The superfluid density as well as the region of superfluid phase diminish as a function of increasing 
temperature. The vanishing of |$o,oP at /2 = and /Z = 10 is an artefact of our approximation 
(see text). 
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FIG. 6: Gritical temperature Tc of the superfluid-normal phase transition as a function of the 
interaction strength U /zt. The solid line is an analytic expression obtained in the approximation 
where we only take into account three slave bosons. The plusses correspond to a numerical solution 
for the case of four slave bosons. 
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FIG. 7: The dispersion relations for k = in the case where we take into account higher filling 
factors at nonzero temperature. On the vertical axis is {hw + ^) / zt and on the horizontal axis is 
JJ . Here we have taken into account all the terms with a = 0, 1, 2 at a temperature of ztl3 = 10. 
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U 

FIG. 8: Dispersion relations ^ + /i as a function of U / zt for k = for zero and nonzero tempera- 
tures. The inner lobe corresponds to zero temperature. The outer lobe correponds to a temperature 
of zt(3 = 3. Here we have only taken into account the first three terms in the right-hand side of 
Eq. ((TH]l. i.e., in the sum we only include the terms with a = and a = 1. 



24 




10 



u 



FIG. 9: The ^-U phase diagram for zero and nonzero temperatures. The inner lobe corresponds 
to the zero-temperature case. The outer lobe corresponds to a temperature of zt^ = 2 




FIG. 10: Total density n as a function of interaction strength U / zt for the first Mott lobe in two 
and three dimensions when including fluctuations. The density approaches a finite value different 
from one, when approaching Uc- 
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